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ABSTRACT 

Motivated by the recent direct detection of cosmological gas infall, we develop an 
analytical model for calculating the mean density profile around an initial overdensity 
that later forms a dark matter halo. We account for the problem of peaks within peaks; 
when considering a halo of a given mass we ensure that this halo is not a part of a 
larger virialized halo. For halos that represent high-sigma fluctuations we recover the 
usual result that such halos preferentially lie within highly overdense regions; in the 
limit of very low-sigma fluctuations, on the other hand, we show that halos tend to lie 
within voids. Combined with spherical collapse, our results yield the typical density 
and velocity profiles of the gas and dark matter that surround virialized halos. 
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1 INTRODUCTION 

Numerical simulations of hierarchical halo forma- 
tion indicate a roughly universal spherically-averaged 
density profile within viriali zed dark matter halos 
iNavarro. Frenk. fc White) Il997f) . Such profiles are gen- 
erated by complicated processes that involve violent 
relaxation and multiple shell crossing; the great interest in 
them results from the possibility of a direct comparison 
to the profiles infe rr ed fro m galactic rotation curves (e.g., 
Sahicci fc BurkertJ 120001 Ide Blok. McGaurii. fc RubirJ 



120011: Ivan den Bosch fc S waters! 1200 lF and from 
lensing in clusters (e.g.. ISand. Treu. fc Ellis! l2002t 
iDahle. Hannestad. fc Sommer-Larserj| ~~ l2002TK However, 
the analogous question of what is the typical density 
profile outside the virial radius has received relatively little 
attention, because of the difficulty of directly observing 
density profiles even out to the halo virial radius; this 
difficulty results from the rapid decrease in the gas and 
dark matter densities with distance from the halo center. 
However, the theoretical description is simpler for infall 
outside the virialization radius, since in cold dark matter 
(hereafter CDM) models there is no shell crossing until the 
final non- linear collapse. 

The gas around halos can only be detected using 
physical effects that are very sensitive to low-density 
gas. One such phenom enon is resonant Lya absorption; 
iGunn fc Peterson! (Il965tl first noted that even an H I neu- 
tral fraction ~ 10" 5 suffices to make gas at the cosmic mean 
density strongly absorb photons at the Lya resonance. Thus, 
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the pattern of gas infall around halos may affect observ- 
able properties of the Lya absorption in these regions. At 
moderate redshifts, these regions near halos play a promi- 
nent role in measurements of the proximity effect of quasar s 
(e.g., IScott et al]l2000t iBaitlik. Duncan, fc Ostrikerlll98gFl . 
In these measurements, the effect of the high quasar inten- 
sity on the statistics of Lya absorption is analyzed and used 
to infer the mean ionizing background, since the quasar in- 
tensity dominates only as long as it is significantly stronger 
than the background. A different application may be pos- 
sible at very high redshifts, corresponding to the time be- 
fore the universe was fully reionized; the presence of the 
fully neutral IGM could be detected through absorption 
of a bright source by the Lya damping wing of the IGM 
jMiralda-Escudel 11998^1 . However, every bright source cre- 
ates a surr ounding H II region which complicated the inter- 
pretation llCen fc Hairnanll2000t livladau fc ReedkOOCTl . The 
properties of the gas in the regions that lie just outside the 
host halo of the source galaxy or quasar should play a major 
role in these various applications of the observed absorption. 

The various calculations of absorption cited above ne- 
glected the enhanced densities and the negative (i.e., in- 
fall) velocities expected for the gas surrounding a massive 
halo at high redshift. This infall, however, may be an es- 
sentia l component of any corre ct interpretation of observa- 
tions; lLoeb fc Eisensteinl il995t) applied it to the proximity 
effect and showed that neglecting infall can lead to an over- 
estimate of the ionizing background flux by up to a factor 
of thr ee. Applying the same infall profile, iBarkana fc Loebl 
(2002) have recently revealed the signature of infalling gas 
around high redshift quasars. The spectral pattern due to 
Lya absorption by the infalling gas can be used to estimate 
the quasar's ionizing intensity and the total mass of its host 
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halo. While the spectra of individual objects are marked by 
the large density fluctuations expected in the intergalactic 
medium (hereafter IGM), if this pattern is averaged over 
many objects then the resulting mean profile could be com- 
pared directly with an accurate theoretical prediction for the 
density and velocity profile of the infalling gas. 



lLoeb fe Eisensteinl (1995) defined a halo in the same 
way as halos were d efined in the classic model of 
iPress fc Schechteil <ll974f) . In this approach, halos on a given 
scale are defined as forming when the initial density aver- 
aged on this scale is higher than a threshold that depends 
on the collapse redshift and is fixed using spherical collapse. 
The halo mass function derived from this assumption, when 
multiplied by an ad-hoc correction factor of two, provides 
a rea sonable match to the result s of numerical simulations 
fe.g.. lKatz. Quinn. fc GeIblll993D . 



The s a me m ass function was later derived by 
iBond et alJ il99ll) using a self-consistent approach, with 
no need for external factors. In their derivation, the factor 
of two has a more satisfactory origin, namely the so-called 
"cloud-in-cloud" problem: if the average initial density on a 
certain scale is above the collapse threshold, this region may 
be contained within a larger region that is itself also denser 
than the collapse threshold; in this case the original region 
should be counted as belonging to the ha lo with mass corre- 
spo nding to the larger collap sed region. IBond et alJ dl99ll) 
and lPress fc Schechterl l)l974h both give the same halo mass 
function and thus appear to be at least statistically valid 
given the agreement with the mass function seen in numer- 
ical simulations. However, the IBond et al.l (|l99lfl model is 
more satisfactory since it accounts for the cloud-in-cloud 
problem and, more importantly, allows for calculations of 
other halo properties in addition to simple mass functions, 
all in an approach that is more self-consistent. 



The IBond et alJ 1^991) approach (also known as 
extended Press-Schechter) h as been used to calcu- 
late halo merger histories llLacev fc Cold I1993T) . halo 
mass functions in models with warm dark matter 
( Bark ana. Haiman. fc Ostrikerl 1200 J ) , halo mass functions 



in a model with ellipsoidal collapse JSheth. Mo. fc Tormerl 
20(n ) , and the nonlinear bias i ng of halos JPorchmT^t^T 
1998; IScannapieco fc Barkanal I2002T) . In this paper, we 



present a new application of this approach in which we ob- 
tain the expected profile of infalling matter around vital- 
ized halos. The rest of this paper is organized a s follows. In 
§21 we establish our notation and review the [ Bond et alJ 
( 1991) derivation of the lPress fc Schechterl B1974T) halo mass 
function. In § 2^2 we r eview the infall profile used by 
1995 



lLoeb fc Eisensteinl J1995J) and also derive an alternative 
model based on similar assumptions but carried out in 
Fourier space. In § 2.3 we derive our model based on the 
extended Press-Schechter formalism. In § 2.4 we then dis- 
cuss the modified picture of halo virialization once infall is 
included. In § 3 we illustrate the predictions of our model 
for the initial and final density profiles surrounding halos. 
Finally in § 4 we give our conclusions. 



2 A MODEL FOR INFALL: DERIVATION 
2.1 Review of Halo Collapse 

Before addressing the problem of the initial density pro- 
file around a virializ ed halo, we briefly review in this sec- 
tion the approach of IBond et alJ il99lT) which leads to the 
halo mass function. We work with the linear overdensity 
field 8(x, z) = /£>(x, z) I p(z) — 1, where x is a comoving po- 
sition in space, z is the cosmological redshift and p is the 
mass density, with p being the cosmic mean density. In the 
linear regime, the density field maintains its shape in co- 
moving coordinates and the overdensity simply grows as 
S = SiD(z)/D(zi), where Zi and Si are the initial reds hift and 
overd ensity, and D(z) is the linear growth factor jPeeblesI 
1980). When the overdensity in a given region becomes non- 
linear, the expansion halts and the region turns around and 
collapses to form a virialized halo. 

The time at which the region virializes can be estimated 
based on the initial linear overdensity, using as a guide the 
collapse of a spherical tophat perturbation. At the moment 
at which a tophat collapses to a point , the overdens ity pre- 
dicted by linear theory is S c — 1.686 llPeebleall980r) in the 
Einstein-de Sitter model. This value depends slightly on cos- 
mological parameters, and even in the Einstein-de Sitter case 
it is modified by infall, as we discuss in H2.4I our analytical 
expressions in this section are valid for any value of S c . 

A useful alternative way to view the evolution of density 
is to consider the linear density field extrapolated to the 
present time, i.e., the initial density field at high redshift 
extrapolated to the present by multiplication by the relative 
growth factor. In this case, the critical threshold for collapse 
at redshift z becomes redshift dependent, 



S c (z) = S c /D(z) 



(1) 



We adopt this view, and throughout this paper the power 
spectrum P(k) refers to the initial power spectrum, linearly 
extrapolated to the present (i.e., not including non- linear 
evolution) . 

At a given z, we consider the smoothed density in a 
region around a fixed point A in space. We begin by aver- 
aging over a large mass scale M, and then lower M until we 
find the highest value for which the averaged overdensity is 
higher than 5 c (z); we assume that the point A belongs to a 
halo with a mass M corresponding to this filter scale. 

In this picture we can derive the mass distribution of 
halos at a redshift z by considering the statistics of the 
smoothed linear density field. If the initial density field is 
a Gaussian random field and the smoothing is done using 
sharp fc-space filters, then the value of the smoothed S un- 
dergoes a random walk as the cutoff value of fc is increased. 
If the random walk first hits the collapse threshold S c (z) at 
fc, then at a redshift z the point A is assumed to belong to 
a halo with a mass corresponding to this value of fc. Instead 
of using fc or the halo mass, we adopt as the independent 
variable the variance at a particular filter scale fc, 

i-k 

(2) 



S k = ^J^ dk'k' 2 P(k') 



In order to construct the number density of halos in this 
approach, we need to solve for the evolution of the probabil- 
ity distribution Q(S, St), where Q(S, Sk) dS is the probability 
for a given random walk to be in the interval 5 to 5 + dS at 
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Sk- Alternatively, Q(S, S%) dS can also be viewed as the tra- 
jectory density, i.e., the fraction of the trajectories that are 
in the interval S to S + dS at Sk , assuming that we consider 
a large ensemble of random walks all of which begin with 

S = at S k = 

iBond et alJ <ll99lT) showed that Q satisfies a diffusion 
equation, 



dQ 
dS k 



l d 2 Q 
2 dS 2 ' 



(3) 



which is solved by the Gaussian solution which we label Qq: 



Qo(5,Sk) = 



V2Ts, 



: exp 



8 2 
~2S k 



(4) 



To determine the probability of halo collapse at a red- 
shift z, we consider the same situation but with an absorb- 
ing barrier at 5 = v, where we set v — S c (z). The fraction 
of trajectories absorbed by the barrier up to Sk corresponds 
to the total mass fraction contained in halos with masses 
higher than the value associated with Sk- In this case, Q 
again satisfies eq. © and the solution with the barrier in 
place is given by adding an extra image-solution: 



Q(u,8,S k ) = Qo(S,Sk)-Qo(2u-S,S k ). 



(5) 



Using this expression, the fraction of all trajectories that 
have passed above the barrier v by Sk is 



F{u,S k 



1 



d5Q{v,5, S k ) 



(6) 



and the differential mass function is found [using eq. to 
be 



f^,S k ) = -^-F(u,S k ) 



2tt5 



3/2 



exp 



v 

2Sk 



(7) 



As f(v, Sk) dSk is the probability that point A is in a halo 
with mass in the range corresponding to Sk to Sk + dSk, the 
halo abundance is then simply 



dn 
dM 



_P_ 

M 



dS k 



dM 



(8) 



where dn is the comoving number density of halos with 
masses in the range M to M + dM. The cumulative mass 
fraction in halos above mass M is similarly determined to 
be 



F(> M\z) = erfc 



(9) 



While these expressions were derived in reference to 
density perturbations smoothed by a sharp fc-space filter as 
given in eq. @, Sk is often replaced in the final results with 
the variance of the mass M enclosed in a spatial sphere of 
comoving radius r: 



a{M) = a(r) = — - / k 2 dkP(k) W 2 (kr) 



2tv 



(10) 



where W(x) is the spherical tophat window function, defined 
in Fourier space as 



W(x) 



sin(:r) cos(x) 



(11) 



and M — 47rp(0)r 3 /3 in terms of the mean density of matter 



at z = 0. With these replacements we rec over the cumulative 
mass fraction that was originally derived iPress fc Schechterl 
Il974l) simply by considering the distribution of overdensities 
at a single point, smoothing with a tophat window function, 
and integrating from S c to 00. In this derivation the au- 
thors were forced to multiply their result by an arbitrary 
factor of two, to account for the mass in underdense re- 
gi ons. The excursion -set derivation presented here, based 
on lBond et alJ il99lT) . properly includes all the mass by ac- 
counting for small regions that lie within overdensities on 
larger scales. This approach also makes explicit the approx- 
imations involved in working with <r 2 (r). Strictly speaking, 
dealing with a real-space filter requires a complete recalcula- 
tion of f(v,a 2 ) which accounts for the correlations intrinsic 
to W(x). However, simply replacing Sk with <r 2 (r) in eq. (|5J 
has been shown to b e in reasonable agreement w ith numeri- 
cal simulations (e.g., Ka tz. Quinn. fc Gelbll993fl . and is thus 
a standard approximation. 

2.2 Infall based on the Press-Schechter Model 

In the Press-Schechter model, a halo of mass M is assumed 
to form at redshift z if the corresponding comoving sphere 
(of ra dius denoted tm) h a s a lin ear overdensity of S = v = 
8c (z). lLoeb fc Eisensteinl jl995l) used this fact to derive a 
mean infall profile as follows. Consider spheres of various 
comoving radii r, all centered on the region in the initial 
density field that contains the mass that ends up in the 
final virialized halo. The mean enclosed overdensity S(r) in 
each such sphere (linearly extrapolated to the present) is a 
Gaussian variable with variance a 2 (r) as given in eq. 1101 . 
The overdensities in two such spheres with radii ri and 
are correlated with a cross-correlation 



&{ri ,r 2 )^ — 



k 2 dkP(k)W(kn)W(kr 2 ) 



(12) 



Thus, the expected value of S(r) given that 8m = S(tm) = v 
is given by their joint Gaussian distribution as 

(5(r))ps- r _ £ r (rM,r) 



o 2 (r M ) 



(13) 



where the notation PS — r indicates a derivation based on 
the Press-Schechter model and carried out in r-space. 

In addition to obtaining the mean profile, it is useful to 
calculate the scatter in the initial profiles around halos of a 
given mass in order to assess the applicability of the mean 
profile to individual objects. If in the future the model is 
to be compared to halos in numerical simulations or in real 
observations, a quantitative formula for the scatter helps 
us determine the number of objects necessary to average 
over before the mean profile emerges. In addition, the scat- 
ter itself is a prediction of the model that in principle can 
be compared to simulations or observations. In the present 
model, the ratio 5(r)/i^ follows a Gaussian distribution with 
the mean value given above and a standard deviation given 
by 



£(r)i 



*(r) &{r M ,r) 



v 2 v 2 a 2 (rM) 



(14) 



Alternatively, we can create an infall model based on 
similar assumptions, but (as in the previous subsection) cal- 
culated in fc-space. Consider a halo that forms on a scale 
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with a corresponding variance Sk,M, and consider a larger 
scale around it with variance Sk < Sk,M- We obtain here a 
Press-Schechter model in fc-space by considering random tra- 
jectories but without including an absorbing barrier; the self- 
consistent inclusion of a barrier, which accounts for peaks- 
within-peaks and is the hallmark of the extended Press- 
Schechter model, is considered in the next subsection. 

In the absence of a barrier, the probability distribution 
of the overdensity S at any Sk is Gaussian, and the values of 8 
at multiple points are jointly normal. Now, 8 has variance Sk 
and Sk.M at the two scales of interest. Furthermore, since the 
random walk from Sk to Sk,M is independent of the random 
walk from to Sk , the cross-correlation of the values of 5 at 
Sk and at Sk,M is given by the overlapping portion of the 
two random walks: the cross-correlation is therefore simply 
Sk- Similarly to the r-space case above, the joint Gaussian 
distribution yields 

(5(r)) PS - k _ S k 



Sk 



(15) 



As in the PS-r model, the distribution of 8(r)/v is Gaussian, 
where now the standard deviation is 



S(r) t 



Sk 



si 



' 2 S k 



(16) 



Note that both of the models developed in this section 
satisfy the continuity condition (8(r))/v — > 1 when r — > tm 
(or equivalently Sk — > Sk,u)- However, we show in the next 
section that they fail to give a physically acceptable result in 
a different limit, when v <§C U Sk.M- Moreover, we have ig- 
nored in this section the failure of the simple Press-Schechter 
approach to account for negative fluctuations, and have not 
included the factor-of-two normalization fudge factor; in the 
next section, in contrast, we develop a model that is self- 
consistent. 



2.3 Infall based on the Extended Press-Schechter 
Model 

In this section we derive a new model for the initial den- 
sity profile around a virialized halo. As additional motiva- 
tion, we first show that the two simple models considered 
in the previous section yield results that are physically un- 
acceptable, once we consider the obvious possibility that a 
small halo will lie within a bigger one. Consider, for exam- 
ple, the PS-r model in a CDM cosmology (the cosmologi- 
cal parameters are listed in section [3J . In order to explore 
the limit of v <C (t(tm) at z = (when v = 1.68 in this 
model), we consider a halo of mass M = 10 6 Mq (tm = 18 
kpc, o{tm) = 8.5), and a surrounding region with r — 2tm 
(M(r) = 8 x 10 6 M Q , a(r) = 7.4). Consider randomly cho- 
sen regions of size r in the universe. Half of them will have a 
negative overdensity and half will have a positive overdensity 
(where, as elsewhere in this paper, we consider the linearly- 
extrapolated overdensity). Those with a positive perturba- 
tion will almost certainly be contained within a halo of mass 
at least M(r), since the typical 8 on this scale is of order 
cr(r), which in turn is much greater than the threshold v 
for forming a halo. Even for those rare perturbations that 
are smaller than v, the fact that the perturbation is positive 
increases the (already very high) chance for the formation 
of some halo with radius between tm and r. Clearly, among 



regions with positive perturbations on the scale r, very few 
will contain an isolated halo of size tm that is not contained 
within a larger halo. 

On the other hand, fully half of such large regions have 
a negative perturbation, again with magnitude of order cr{f). 
Such a negative perturbation implies that certainly no halo 
has formed on the scale r; together with the correlation 
among spheres of different radii, it also implies a very strong 
reduction in the chance that any halo has formed on a scale 
between tm and r. This double effect implies that isolated 
halos of size tm will be found almost exclusively within voids 
on the scale r. Indeed, in the limit of v /cr(rjvf) — > 0, 8 on the 
scale r should be negative with a magnitude at least com- 
parable to er(r/vf), in order for 8 to remain negative at all 
scales from r down to tm- Thus, in this limit we should find 
that the ratio S(r)/u — > — oo. This contrasts with the Press- 
Schechter models from the previous section, which predict 
that the typical such case should have a positive 8 on the 
scale r (in fact, a 8 only slightly below i>); this behavior 
arises from the general failure of these models to properly 
incorporate negative density fluctuations. 

To derive our new, self-consistent model, we apply the 
extended Press-Schechter model in Fourier space and then 
use an ansatz to convert it to real space. We begin by follow- 
ing the derivation of (<5(r))ps_k in the previous subsection 
but now including the absorption barrier at 8 = v. Thus, 
we consider once again the two scales with variances Sk.M 
and Sk, but we now require Sk,M to form a halo that is not 
contained within any larger virialized halo; in particular this 
implies that there is no halo on the scale corresponding to 
Sk or on any larger scale. Then the probability distribution 
of the overdensity 8 at Sk is given by Q(y, 5, Sk) as in eq. JjJ. 
As before, the random walk from Sk onwards is independent 
of the random walk leading up to Sk- In order to account 
for the extended Press-Schechter definition of halo forma- 
tion ; we use argument s sim ilar to those leading to eq. Q 
[see lBond et alJ Jl99 Jl and lLacev fc Cold <1993ft for other 
applications of similar arguments]; the probability distribu- 
tion of the initial density at Sk, given the formation of a 
halo at Sk ai, is found to be 



P(6\6i 



Q(v,6,Sk) f{l/ -j ,S o M : Sk) ■ (17) 



f{y,Sk, M ) 

The mean expected 8 is given by integrating <5 over this 
distribution, with limits of — oo to +v; we thus obtain the 
result of using the extended Press-Schechter model in k- 
space: 



-a + 



erf 



13(1- 



2a(l 



7T/J 



exp 



2a 



/J(l -a) 



2a 



where we have defined 



Sk 

Sk.M 



Sk. 



(18) 



(19) 



Unlike the models in the previous subsection, the 
present model does not yield a Gaussian probability distri- 
bution for the density. To find the scatter, we first integrate 
the probability distribution given by eq. 117H and find the 
cumulative probability distribution Pc- Using the variable 
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7 = 5/^ instead of 5, the probability of obtaining a value of 
8(r) that is between — oo and vy equals 



(5(r)) _ (5(r)) cPS _ 



(22) 



1 + - [erf (A) + eri(B)} 



+ 



271-/3(1 - a) 



where we use the shorthand 



A 



( 7 + a-2)V? 



(20) 



(21) 



We can now numerically calculate the 1—a interval of 5(r) by 
bracketing the central 68.3% of the probability. For a given 
M and v, we thus define the minus 1 — a density profile as 
given by values of 7 determined at each r by Pc(t) = 15.9%, 
while the plus 1 — a profile is determined by Pail) = 84.1%. 

Having developed several possible infall models, we 
compare their relative merits by considering several limit- 
ing cases. First, when r — > rjw (or equivalently 5fc — > Sk,i\i), 
all the models satisfy the continuity condition (5{r))/v — > 
1 with zero scatter. In the limit of /3 — *• 0, we require 
{S(r))/v — > —00 as explained earlier in this section. While 
(<5(r)) ps-k behaves correctly in this limit (and approaches 
a constant that does not depend on v), (5(r))ps_ r and 
(<5(r))ps_k are proportional to v and thus fail to correctly 
describe the j3 limit. Thus, only (5(r)) e ps-k behaves 
physically in both of the limits r — > rM and j3 0. 

In the limit of a halo that represents an extremely rare 
peak, i.e., j3 — > 00, the Press-Schechter model becomes an ac- 
curate description that is equivalent to the extended Press- 
Schechter model; the reason is that if the central halo is 
extremely rare then it is exceedingly improbable to find a 
virialized halo on any scale larger than rjw, and thus it is 
acceptable to neglect the barrier at <5 = v except right near 
the scale tm- In this limit, (<5(r)) e ps-k and (<5(r , ))ps-k both 
approach the value va. However, in this limit only, the PS-r 
model is the best available model since it directly considers 
spheres in real space, and these correspond most directly to 
the regions that actually collapse to form halos; this model 
also applies most directly to the threshold v that is derived 
based on spherical collapse in real space. Indeed, the best 
possible model would be (5(r)) e ps-r, a quantity based on 
the extended Press-Schechter model applied directly in real 
space. However, such a model would be extremely difficult to 
develop, since excursion sets of 5 as a function of a 2 (r) would 
no lon ger correspond to pure random walks [see lBond et alJ 
lll99lf) for related discussion]. Thus, following the tradition 
of the extended Press-Schechter approach (see £ffl and l2"Tt . 
we use (<5(r)) e ps-k along with an ansatz. 

To arrive at the most logical ansatz, we simply replace 
a and f3 in eq. 11911 with their r-space equivalents. Based 
on the above arguments, we adopt (<5(r))ps- r and E(r)ps- r 
as the correct j3 — > 00 limits of the mean and the standard 
deviation, respectively. Indeed, 5(r) in the ePS-k model ap- 
proaches a Gaussian distribution in the limit of (3 — > 00, 
and therefore we can match this limit identically to the PS-r 
model, and this procedure fixes unambiguously the required 
r-space versions of both a and j3. Thus, our final ansats for 
the mean initial profile of the density averaged in spheres of 
radius r, linearly extrapolated to z — 0, around a virialized 
halo is 



v v 
where we use eq. (1181 except that we define 

a _ &{r M ,r) _ _ v 2 a{l-a) 

a 2 (r M ) ' a 2 (r) - a£ r (r M ,r) 

We apply these same definitions to eqs. 1201 and 12 1 1 in 
order to estimate the scatter in the infall profiles as given 
by this final model. 

2.4 Halo Virialization in the Presence of Infall 

As noted in H2.ll the linearly extrapolated initial overden- 
sity of a sp herical tophat corresponding to the moment of 
collapse is JPeeblesHl98ol) S c = 3(12tt) 2/3 /20 if Q, m = 1. In 
reality, even a slight violation of the exact symmetry of the 
initial perturbation can prevent the tophat from collapsing 
to a point. Instead, the halo reaches a state of virial equilib- 
rium by violent relaxation (phase mixing). Using the virial 
theorem U = —2K to relate the potential energy U to the 
kinetic energy K in the final state, the final overdensity rel- 
ative to the co smic mean density at the collapse redshift is 
JPeeblesHl980l) A c = 18tt 2 in the Einstein-de Sitter model. 
Thus, if a tophat containing mass M is predicted by spheri- 
cal collapse to reach a physical radius R = at the collapse 
redshift 2, this amount of mass is instead assumed to be dis- 
tributed within a sphere of size equal to the (physical) virial 
radius 7i v i r determined by M = ^■Kp(z)A c R^ il . 

In the presence of infall, however, the enclosed mass is 
modified. Indeed, at the moment when the tophat mass M 
is predicted by spherical collapse to reach a radius R = 0, 
the surrounding infall implies that additional mass has al- 
ready accumulated within the sphere of radius i? v ir- Thus, 
the overdensity within the standard virial radius i? v j r is pre- 
dicted to be higher than the standard value of 187r 2 . There- 
fore, once infall is included, even in the Einstein-de Sitter 
model the precise overdensity to be used to define a "viri- 
alized" halo becomes somewhat uncertain. In this paper we 
choose the overdensity of 18-7T 2 to define the virial radius 
in all cosmological models, since this overdensity has been 
frequently used to compare the analytical model predictions 
(e.g., for the halo mass fu nction) with numerical simulations 
(e.g.. IJenkins et alJl2001f) . Thus, we carry out spherical col- 
lapse calculations with the average infall profile derived in 
the previous section. Having fixed the virial overdensity we 
then find the enclosed virial mass M v j r , which is greater 
than the original tophat mass M by up to a factor of two 
(depending on M, z, and the cosmological parameters). We 
also find the linearly extrapolated initial overdensity of the 
virialized mass shell, which we denote S v . 

In the Einstein-de Sitter case we find S v = 1.5927 for 
the infall-based definition of virialization just noted, com- 
pared to the standard value of 5 C = 1.6865 for an iso- 
lated initial tophat perturbation. More generally, 5 C de- 
pends on the cosmological parameters. In a flat cosmology 
with Q m (z) + £Ia{z) = 1 at a redshift 2, we find the fol- 
lowing fitting function accurate to 5 x 10~ 5 in the ran ge 
0.05 < fi m (2) < 1 (see also, e.g.. lEke. Cole fc Frenklll996h : 



1.6865 



1 + 



1000 



(4.86 - .51a; - .113x 



(24) 



where x = log [fl m (.z)]. Although S v does not strictly depend 
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only on the cosmological parameters, we find numerically 
that it is essentially independent of halo mass. Indeed, over 
a wide parameter range of 10 4 A/q < M < 1O 16 M0 and 
0.05 < Sl m (z) < 1, S v is given to within a relative accuracy 
of 3 x 10~ 4 by 



S v = 1.5927 



1 



1000 



3.9 - .36x - .09x 



(25) 



Thus, for example, 8 V = 1.573 for il m (z) = 0.05, compared 
to 5 C = 1.659 in this case. 

We note that the Press-Schechter halo mass function is 
known to underestimate the abundance of rar e halos relative 
to th e results of numerical simulations (e.g., I.Tenkins et alJ 
l200lT) . This can be fixed in the Einstein-de Sitter model if 
the collapse threshold is multiplied by a factor of 0.87, while 
we find that infall lowers the threshold, from 5 C to S v , by a 
factor of only 0.94. Thus, infall by itself cannot fully explain 
the discrepancy in the Press-Schechter prediction for the 
number of high-mass halos. 



3 RESULTS 

We illustrate results for the standard cosmological param- 
eters Q m = 0.3, Qa = 0.7, h = 0.7 and a CDM power 
spectrum (with initial slope n = 1) normalized to a variance 
as — 0.8 in spheres of radius 8 /i -1 Mpc. We consider a 
wide range of redshifts but focus on the higher values since 
at z < 2 the infalling intergalactic gas is very tenuous and 
may thus be undetectable; however, note that at low redshift 
a significant fraction of the infalling baryons may already be 
in the form of virialized dwarf galaxies that are themselves 
detectable. 

Figure Q shows the predictions of our model [eqs. I|22|l 
and (23)] for the mean initial profiles of enclosed density 
around virialized halos. In most cases, the profiles lie signif- 
icantly below the asymptotic case of a rare halo (/3 — > oo) 
for which the profile is given by (5(r))ps- r as in eq. (13) ■ In 
particular, when M v i T = 1O 8 M0, a halo forming at z < 2 is 
most likely to be found in a void on all scales larger than a 
few times tm- Note that we restrict all calculations to radii 
where the initial profile of enclosed density is monotonically 
decreasing, since otherwise shells would cross during the sub- 
sequent collapse; this cutoff affects only the z = cases as 
well as the M v i r = 1O 8 M0 halo forming at z — 2 or 3. 

We calculate spherical collapse with the initial profiles 
shown in Figure Q The resulting final density profiles are 
shown in Figure [5] In the limiting case of an extremely 
rare halo (i.e., very high formation redshift at a fixed halo 
mass), low mass halos have higher densities of dark mat- 
ter surrounding their virial radius. This can be understood 
roughly as follows. For a power spectrum of density fluctu- 
ations approximated over some range of scales as a power 
law P(k) oc k n , the correlation function over the same range 



of scales varies as f (r 



-(n+3) 



The initial density profile 



around a rare halo follows, roughly, the dark matter corre- 
lation function. On small scales, the CDM power spectrum 
approaches n = — 3 and the initial density profile is rather 
flat, while on larger scales the power index increases toward 
n — —2 for large galaxies, which are therefore surrounded by 
an initial profile that falls off more steeply (see Figure 0. 
These differences are magnified by the process of gravita- 
tional collapse, resulting in final values for p/p at the virial 




-0.4 



r/r M 



10 



30 



Figure 1. Mean initial density profiles (linearly extrapolated to 
the present) around virialized halos. Note that <5(r) is averaged 
over a sphere of comoving radius r, and rjvf is the initial comoving 
radius of the tophat. Results are shown for final halo virial masses 
M vir = 1O 8 M (bottom panel), 10 10 Af Q (middle panel), and 
10 12 Mq (top panel). In each case, profiles are shown (solid curves) 
at z = 0, 2, 3, 5, and 10, respectively, from bottom to top. Also 
shown is the asymptotic case of an extremely rare halo (dashed 
curve). The position of zero overdensity is also indicated (dotted 
line). 



radius of 41, 34, and 27 for M viT = 10 8 , 10 10 , and 1O 12 M , 
respectively (see Figure |2). However, at moderate redshifts 
high mass halos represent rarer, more extreme peaks than 
do low mass halos. Thus, e.g., even at redshift 5 the profile 
surrounding M v i r < 1O 1O M0 halos falls significantly below 
the profile of the extreme limiting case. 

Next, we use eq. 12011 to predict the expected scatter in 
the density profiles around virialized halos. Figure [3] shows 
the initial and final density profiles for M v j r = 1O 1O M0 at 
two different redshifts, for the mean expected profile that we 
have thus far considered as well as for the plus and for the 
minus 1 — a profiles. The final density at J? v ir, for example, 
is higher than the mean density by a factor of 12j^| at z = 2 



and 28" 1 



at z 



10. We thus find a significant scatter, 



although it is smaller at the higher redshift. If the scatter 
as calculated within our spherical model is representative of 
the true scatter among halos, its magnitude suggests that 
our predicted mean profile can be compared accurately to 
an observed sample if the sample contains at least a few 
dozen halos of comparable masses and redshifts. 

Figure0]shows the final profiles of peculiar velocity V pe c 
of the infalling dark matter corresponding to the halos con- 
sidered in Figures Q and |5] The same trends are seen in 
this figure as in the previous ones; in the extreme case of a 
rare halo, infall is stronger for low halo masses, but at mod- 
erate redshifts, the higher-mass halos actually produce the 
stronger infall. Figure0]also shows the peculiar velocity that 
would correspond to stationary matter, i.e., with zero total 
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Figure 2. Final density profiles around virialized halos with the 
initial profiles shown in Figure Note that p(R) is the actual 
density (not enclosed average) at a physical radius R, p is the 
cosmic mean density at the virialization redshift, and i? v j r is the 
final virial radius defined as enclosing a mean density of 187r 2 
times p. Results are shown for final halo virial masses M v j r = 
10 s M (bottom panel), 10 10 M Q (middle panel), and 10 12 M o 
(top panel). In each case, profiles are shown (solid curves) at 
z = 0, 2, 3, 5, and 10, respectively, from bottom to top. Also 
shown is the asymptotic case of an extremely rare halo (dashed 
curve). The position of the cosmic mean density is also indicated 
(dotted line). 

velocity (long-dashed curve). This curve is given by minus 
the Hubble expansion velocity VHub, which in these units is 

V^IR 

Vc J7T fivir 

Values of V pcc that lie below — Vkub give a total velocity 
corresponding to infall. 

When gas falls into a dark matter halo, incoming 
streams from all directions strike each other at supersonic 
speeds, creating a strong shock wave. Three-dimensional hy- 
drodynamic simulations show that the most massive ha- 
los at any time in the universe are ind eed surrounded by 
strong, quasi-spherical accreti on shocks l|Keshet et al]l2o"o2t 
I Abel. Brvan. fc Normanll2002l) . The radius of this shock de- 
pends in general on the history of gas infall into the halo, 
and will thus vary somewhat among halos depending on the 
history of their assembly through mergers. In particular, in 
the spherically-symmetric secondary infall solution, when a 
trace amount of adiabatic gas is allowed to flow in a grav- 
itational potential determined by collisionless dark matter 
with density f2 m = 1, the accretion shock forms at the ra- 
dius -R90.8 which encloses a mean dark matter density of 
90.8 times the c osmic mean <|B ertschin ge r l| T985). However, 
the simulation of lAbel. Brvan. fc Normanl j2002( ) shows that 
along most lines of sight the shock occurs at a smaller dis- 
tance, close to the virial radius. 

The position of the accretion shock is an important 



Figure 3. Scatter in the density profiles around virialized halos. 
Shown are the initial enclosed density profile (top panel) and 
the final density profile (bottom panel), defined as in Figures HI 
and |21 respectively. Each curve fixes the final halo virial mass to 
Af v ; r = 10 10 Mq, with the halo forming at z = 2 (solid curves) or 
z = 10 (dashed curves). In each case, profiles are shown for the 
minus 1 — a profile, the mean expected profile, and the plus 1 — a 
profile, respectively, from bottom to top (see text). 



theoretical prediction that can be probed by observations, 
since the infalling gas that is about to cross the shock front 
produce s a sharp absorption fea ture in high-redshift quasar 
spectra fearkana fc Loebl 120021) . In order to illustrate the 
range of possible properties of this gas element that lies 
just beyond the shock front, we use the results mentioned 
above as a guide and consider two cases for the shock ra- 
dius Rsh'- the smaller radius R s h = R v i r or the larger radius 
R B h — Rgo.s- We assume that the pressure gradient is negli- 
gible in the pre-shock gas compared to gravity, and therefore 
the gas density and velocity are equal to those of the dark 
matter at distances beyond the shock front. 

Figure |5] shows the density and total infall velocity of 
the gas that lies just beyond the shock front. The density 
and the magnitude of the infall velocity both tend to in- 
crease with redshift, and follow only weak trends with halo 
mass at a given redshift. The differences between the two 
redshifts shown are fairly significant for low mass halos but 
rather small for high mass halos. In general, the infalling 
gas is significantly overdense even for halo masses corre- 
sponding to z = dwarf galaxies, especially at the high 
redshifts which are increasingly becoming the focus of ob- 
servations. The figure allows us to compare two sources of 
scatter, namely the varying position of R s h and the variabil- 
ity in the initial density profile. The two effects produce a 
comparable scatter in the pre-shock density, while the in- 
fall velocity, which depends on the enclosed mean density 
rather than on the density itself, is primarily sensitive to 
the position of the shock. 
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R/R»lr 

Figure 4. Final velocity profiles around virialized halos with the 
initial profiles shown in Figure^ Here V pec is the peculiar velocity 
(where negative values correspond to infall) at a physical radius 
R, and V c is the circular velocity at the halo virial radius i? v j r . 
Results arc shown for final halo virial masses M v ; r = 10 s Mq (bot- 
tom panel), 1O 1O M0 (middle panel), and 1O 12 M0 (top panel). In 
each case, profiles are shown (solid curves) at z = 0, 2, 3, 5, and 
10, respectively, from top to bottom. Also shown is the asymp- 
totic case of an extremely rare halo (short-dashed curve). For 
comparison, the curve corresponding to uniform cosmic expan- 
sion is shown (dotted line), as well as the curve corresponding to 
stationary matter (long-dashed curve). 
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Figure 5. Properties of the infalling gas that is just about to 
cross the accretion shock around a virialized halo. Shown are the 
density in units of the cosmic mean (top panels), and the to- 
tal (not peculiar) velocity in units of the halo circular velocity 
(bottom panels), both as a function of the final halo virial mass 
Af v i r . Results are shown at z = 2 (left panels) and z = 10 (right 
panels). Two possible values are considered for the shock radius, 
R B h = -Rvir (dashed curves) and R s ^ = -R90.8 (solid curves). In 
each case, profiles are shown for the minus 1 — a profile, the mean 
expected profile, and the plus 1 — a profile, respectively, from bot- 
tom to top in the top panels (the order is reversed in the bottom 
panels). 



4 CONCLUSIONS 

We have developed a model for calculating the initial den- 
sity profile around overdensities that later collapse to form 
virialized halos. This is the first such model to account for 
the possibility that an overdensity on a given scale may be 
contained inside a large overdensity on an even bigger scale. 
As a result, we have found that the mean expected profile 
[eqs. 118H and 12311 1 depends on both the mass of the halo 
and its formation redshift. Starting from this mean initial 
profile, we have used spherical collapse to obtain the final 
density profile surrounding the virialized halo. 

In reality, there will be some variation in the initial 
profiles, with each leading to a different final profile. We 
have estimated the scatter within our spherical model, and 
found that a sample of around a few dozen halos is required 
in order to obtain an accurate mean profile by averaging. 
Halo samples derived from numerical simulations or from 
observations can be used to test our predictions for the mean 
and for the scatter of the density and velocity profiles of 
infalling matter around virialized halos. 

Once infall is included, the standard derivation of the 
mean density enclosed within the virial radius fails. Redefin- 
ing the virial radius as the radius enclosing a density higher 
than the cosmic mean by the standard value of 187T 2 , we 
have found that the initial overdensity needed for a halo to 
form at some final redshift is lower than it would have to be 



for an isolated tophat perturbation [compare eqs. 1241 and 
125l ll . Therefore, accounting for infall increases the predicted 
abundance of rare halos for a given initial power spectrum, 
but the increase is too small to fully explain the discrepancy 
between the Press-Schechter prediction for the number of 
high-mass halos and the number seen in numerical simula- 
tions. 

We have confirmed that rare halos at a given redshift 
are surrounded by large, overdense regions of infall, but we 
have found that the more numerous halos that correspond to 
low-sigma peaks are surrounded by small infall regions that 
lie within large voids. However, the infalling gas just beyond 
the cosmological accretion shock is in general much denser 
than the cosmic mean; a density larger by a factor of at least 
10 occurs around the most massive halos at z > 2. Possible 
applications of our results include the Lyman-a absorption 
profiles of high-redshift objects and the proximity effect of 
quasars at all redshifts. 
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